Scale-Dependent  Infrared  Radiative  Damping  Rates  on  Mars  and 
Their  Role  in  the  Deposition  of  Gravity -Wave  Momentum  Flux 

Stephen  D.  Eckermann^’*,  Jun  Ma’’,  Xun  Zhu‘^ 

“Space  Science  Division,  Code  7646,  Naval  Research  Laboratory,  4555  Overlook  Avenue  SW,  Washington,  DC 

20375,  USA 

^Computational Physics,  Inc.,  8001  Braddock  Road,  Springfield,  Virginia  22151,  USA 
“Johns  Hopkins  University/Applied  Physics  Laboratory,  11100  Johns  Hopkins  Road,  Laurel,  Maryland  27023,  USA 


Abstract 

Using  a  Curtis  matrix  model  of  15  yum  CO2  radiative  cooling  rates  for  the  Martian  atmosphere,  we 
have  computed  vertical  scale-dependent  IR  radiative  damping  rates  from  0-200  km  altitude  over 
a  broad  band  of  vertical  wavenumbers  m  =  2;r(  1-500  km)"'  for  representative  meteorological 
conditions  at  40°N  and  average  levels  of  solar  activity  and  dust  loading.  In  the  middle  atmo¬ 
sphere,  infrared  (IR)  radiative  damping  rates  increase  with  decreasing  vertical  scale  and  peak  in 
excess  of  30  days"'  at  ~50-80  km  altitude,  before  gradually  transitioning  to  scale-independent 
rates  above  ~100  km  due  to  breakdown  of  local  thermodynamic  equilibrium.  We  incorporate 
these  computed  IR  radiative  damping  rates  into  a  linear  anelastic  gravity-wave  model  to  assess  the 
impact  of  IR  radiative  damping,  relative  to  wave  breaking  and  molecular  viscosity,  in  the  dissipa¬ 
tion  of  gravity-wave  momentum  flux.  The  model  results  indicate  that  IR  radiative  damping  is  the 
dominant  process  in  dissipating  gravity-wave  momentum  fluxes  at  ~0-50  km  altitude,  and  is  the 
dominant  process  at  all  altitudes  for  gravity  waves  with  vertical  wavelengths  <10-15  km.  Wave 
breaking  becomes  dominant  at  higher  altitudes  only  for  “fast”  waves  of  short  horizontal  and  long 
vertical  wavelengths.  Molecular  viscosity  plays  a  negligible  role  in  overall  momentum-flux  depo¬ 
sition.  Our  results  provide  compelling  evidence  that  IR  radiative  damping  is  a  major,  and  often 
dominant  physical  process  controlling  the  dissipation  of  gravity  wave  momentum  fluxes  on  Mars, 
and  therefore  should  be  incorporated  into  future  parameterizations  of  gravity-wave  drag  within 
Mars  GCMs.  Lookup  tables  for  doing  so,  based  on  the  current  computations,  are  provided. 
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1.  Introduction 


Being  ~95%  CO2  by  mixing  ratio,  the  thermal  balance  of  the  Martian  atmosphere  is  driven 
to  a  large  degree  by  infrared  (IR)  radiative  transfer  (RT)  among  CO2  bands,  particularly  15  //m 
cooling  and  heating  by  absorption  of  solar  radiation  in  near-IR  bands.  The  advent  of  satellite 
remote  sensing  of  the  Martian  atmosphere  has  spurred  development  of  detailed  atmospheric  RT 
codes  (Lopez- Valverde  and  Lopez-Puertas,  1994),  to  improve  both  the  retrieval  of  atmospheric 
parameters  from  passive  IR  remote  sensors  (Kleinbohl  et  ah,  2009),  as  well  as  the  parameteriza- 
tions  of  radiative  heating  and  cooling  rates  embedded  within  Martian  general  circulation  models 
(GCMs:  e.g.,  Gonzalez- Galindo  et  ah,  2009).  These  models  have  highlighted  the  importance  of 
departures  from  local  thermodynamic  equilibrium  (LTE)  in  CO2  IR  RT  that  significantly  modify 
net  heating  and  cooling  rates  above  ~80  km  altitude  (Lopez-Puertas  and  Lopez- Valverde,  1995; 
Lopez- Valverde  et  ah,  1998). 

Finite  computational  resources  make  it  impractical  to  incorporate  detailed  RT  codes  within 
Martian  GCMs  to  compute  net  CO2  IR  heating  and  cooling  rates,  and  much  simpler  (i.e.  faster) 
parameterizations  are  needed.  CO2  IR  radiative  cooling  rates  are  currently  parameterized  in  Mars 
GCMs  at  high  altitudes  by  local  cooling-to-space  approximations  tuned  to  reproduce  results  of 
detailed  offline  non-LTE  RT  calculations,  and  at  low  altitudes  by  simplified  wide-band  ETE-based 
schemes  (Forget  et  ah,  1999;  Richardson  et  ah,  2007;  Gonzalez-Galindo  et  ah,  2009). 

These  same  computational  constraints  force  the  GCMs  to  operate  at  gridpoint  resolutions  that 
are  too  coarse  to  resolve  the  full  spatial  spectrum  of  atmospheric  dynamics  that  can  exist  in  the 
Martian  atmosphere  (Imamura  et  ah,  2007).  Gravity  waves,  for  example,  are  mostly  unresolved 
by  Martian  GCMs,  and  so  the  missing  effects  of  their  momentum  flux  deposition  on  the  resolved 
circulation  (so-called  gravity-wave  drag)  must  be  parameterized,  again  using  simplified  fast  algo¬ 
rithms  (e.g.,  Joshi  et  ah,  1996;  Collins  et  ah,  1997;  Forget  et  ah,  1999). 

Focal  IR  cooling  parameterizations  ignore  the  inherently  nonlocal  nature  of  radiative  cooling: 
for  example,  that  an  atmospheric  layer  wedged  between  two  colder  layers  will  generally  cool  more 
rapidly  than  if  the  temperature  is  identical  in  all  three  layers.  Thus  vertical  structure  in  atmospheric 
temperature  profiles  can  significantly  change  IR  cooling  rates.  At  higher  altitudes  where  cooling 
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rates  on  Mars  are  the  largest,  most  of  the  short  vertical  scale  structure  in  atmospheric  temperature 
profiles  is  produced  by  gravity  wave  oscillations  (Fritts  et  ah,  2006;  Creasey  et  ah,  2006),  which 
the  GCMs  do  not  resolve.  This  suggests  a  potentially  important  interaction  between  gravity  wave 
dynamics  and  IR  radiative  cooling  that  the  existing  uncoupled  parameterizations  of  each  process 
in  Mars  GCMs  miss  entirely. 

The  general  issue  has  been  recognized  in  planetary  atmosphere  for  many  years  (e.g.,  Spiegel, 
1957;  Sasamori  and  London,  1966;  Goody  and  Belton,  1967).  In  the  Earth’s  atmosphere.  Eels 
(1982)  showed  that  standard  WKB  methods  for  gravity  waves  could  be  applied  in  the  computation 
of  the  gravity-wave-modulated  IR  cooling  rate,  which  yielded  a  local  IR  radiative  relaxation  rate 
that  acted  to  damp  the  gravity  wave-induced  temperature  perturbations.  This  damping  rate  varied 
according  to  the  wave’s  vertical  wavelength  A^.  The  Eels  approach  has  been  extended  in  various 
ways  to  more  realistic  non-ETE  RT  models  (Shved  and  Utyakovsky,  1983;  Eels,  1984;  Zhu  and 
Strobel,  1991;  Zhu,  1993;  Bresser  et  ah,  1995)  and  the  resulting  scale-dependent  IR  damping  rates 
have  been  incorporated  within  parameterizations  of  gravity  wave  dynamics  on  Earth  (e.g.,  Marks 
and  Eckermann,  1995;  Eckermann  and  Marks,  1997). 

The  analytical  methods  of  Eels  (1984)  were  applied  by  Imamura  and  Ogawa  (1995)  to  15  /rm 
CO2  cooling  within  the  atmospheres  of  Earth,  Mars  and  Venus.  On  Mars,  they  inferred  scale- 
dependent  IR  damping  rates  that  were  significantly  faster  than  on  Earth.  They  also  estimated  the 
relative  impacts  of  molecular  viscosity  and  radiative  cooling  on  the  dissipation  of  gravity  waves  as 
a  function  of  altitude  and  A^  on  each  planet.  Despite  these  suggestive  preliminary  findings,  there 
has  been,  as  far  as  we  are  aware,  no  additional  research  on  scale-dependent  IR  radiative  damping 
on  Mars.  In  particular,  it  continues  to  be  unclear  whether  IR  radiative  damping  has  any  relevance 
to  the  dissipation  of  gravity-wave  momentum  flux  on  Mars. 

Thus,  this  study  has  two  parts.  Eirst,  in  section  2  we  formulate  and  validate  a  Curtis-matrix 
model  of  15  //m  radiative  cooling  on  Mars,  and  apply  it  in  section  3  to  compute  scale-dependent 
IR  radiative  damping  rates.  Then,  in  section  4  we  incorporate  these  damping  rates  into  sim¬ 
ple  anelastic  models  of  gravity  wave  dynamics  on  Mars,  to  assess  their  impact  on  gravity-wave 
momentum-flux  deposition,  relative  to  wave  breaking  and  molecular  viscous  dissipation.  Section  5 
summarizes  our  major  findings  and  conclusions  concerning  the  impact  of  IR  radiative  damping  on 
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gravity- wave  momentum-flux  deposition  on  Mars. 


2.  CO2  15  /im  Cooling  Rates  for  Mars 

2.1.  Model 

Exaet  ealeulations  of  CO2  IR  eooling  rates  involve  solving  the  RT  equation  numerieally  using 
first-prineiples  line-by-line  integration  through  an  optieally  thiek  Martian  atmosphere.  High  ae- 
euraey  requires  multiple  integrations  over  a  range  of  radianee  paths  and  a  fine  IR  frequeney  grid, 
whieh  is  prohibitively  expensive  eomputationally,  and  thus  performed  only  oeeasionally  for  vali¬ 
dation  purposes  (e.g.,  Zhu,  1990;  Fomiehev,  2009).  In  eommon  with  many  previous  studies  of  IR 
eooling  on  Earth  (e.g.,  Zhu,  1990)  and  Mars  (e.g.,  Eopez-Valverde  and  Eopez-Puertas,  1994),  we 
eompute  a  faster  approximate  solution  using  Curtis  matriees. 

Ignoring  boundary  flux  terms,  the  eooling  rate  expressed  in  Curtis-matrix  form  is 

dT 

Q  ^  —  =  C0.  (1) 

dt 

Here  Q,  T  and  0  are  eolumn  veetors  representing  a  diseretized  vertieal  profile  of  eooling  rate, 
temperature  and  normalized  souree  funetion  at  a  series  of  vertieal  grid  points  k  =  1 ...  L  at  sueees- 
sive  altitude  (pressure)  levels  zu  (Pk),  and  C  is  the  L  x  L  Curtis  matrix  that  speeifies  the  generally 
nonloeal  dependenee  of  eooling  rate  on  temperature  and  mixing  ratio. 

Inelusion  of  a  eomplete  list  of  CO2  vibration-translation  (V-T)  and  vibration- vibration  (V-V) 
transitions  in  (1)  yields  a  large  number  of  nonlinearly  eoupled  rate  and  RT  equations  whieh  ean 
be  very  diffieult  or  even  impraetieal  to  solve  numerieally  (Kutepov  et  ah,  1998).  For  the  redueed 
eolleetion  of  major  V-T  and  V-V  transitions  used  by  Dickinson  (1984),  Zhu  (1990)  showed  that 
nonlinearity  in  the  resulting  equations  arose  from  coupling  terms  due  to  V-V  exchanges  which,  if 
eliminated,  reduced  the  system  to  an  equivalent  two-energy-level  problem.  While  V-V  transitions 
are  clearly  significant  to  the  RT  in  individual  sub-bands,  Zhu  (1990)  argued  that  they  were  less 
significant  in  determining  total  non-ETE  cooling  rates  over  the  entire  band.  Offline  calculations 
by  Zhu  and  Strobel  (1990)  found  errors  of  no  more  than  ~1  K  day^^  in  total  IR  cooling  rate  at 
the  Earth’s  mesopause  when  V-V  transitions  were  ignored.  For  Mars,  Eopez-Puertas  and  Eopez- 

Valverde  (1995)  also  found  little  change  in  total  cooling  rate  when  V-V  rates  were  varied.  Their 
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results  also  indicated  that  the  contribution  of  V- V  transitions  to  the  total  cooling  rate  in  the  Martian 
atmosphere  is  smaller  than  on  Earth.  Hence  we  ignore  V-V  transitions  in  what  follows. 

For  the  resulting  equivalent  two-energy-level  Curtis-matrix  model  of  Zhu  (1990),  Zhu  et  al. 
(1992)  developed  an  algorithm  for  deriving  cooling  rates  for  arbitrary  atmospheres  by  interpolating 
a  reference  Curtis  matrix  Co,  computed  for  reference  profiles  of  pressure  po,  temperature  Tq  and 
CO2  mixing  ratio  to  yield  a  solution  for  given  p,  T  and  The  element  {i,k)  of  the 

interpolated  Curtis  matrix  C  is  derived  by  modifying  the  reference  value  as 

C,.  =  +  A'«ST,  + 

+  A‘i‘(<5r,)(v„)  +  y«ln„  +  (2) 


The  matrices  defined  in  (24a)-(24e)  of  Zhu  (1990),  are  differencing  approximations  of  deriva¬ 
tives  resulting  from  Taylor  series  approximation  of  the  escape  function  at  levels  pj  and  pk,  and  6Tj 
is  a  temperature  correction  from  the  reference  (Tq)  to  actual  temperature  T  at  level  pj.  Following 
Fels  and  Schwarzkopf  (1981),  the  matrix 


Y}n=jGn5Tn 
YL=j  Gn 


(3) 


in  (2)  encapsulates  the  change  in  transmission  due  to  an  integrated  change  in  temperature  between 
levels  Pj  and  p^,  based  on  the  precomputed  weight  vector  G  for  the  reference  state.  The  last  two 
terms  in  (2)  with  coefficient  matrices  and  y^j  represent  the  first-  and  second-order  corrections 
due  to  the  effects  of  integrated  changes  in  CO2  column  density  between  levels  j  and  k,  where  (Zhu, 
1994) 

l>Qn  ) 

ii/T  —  ,  (4) 

ZLyG„(ln2-ln0.5)/2 

Non-FTE  correction  of  the  Curtis  matrix  (2)  follows  the  equivalent  two-level  approach  of  Zhu 
(1990)  as  summarized  by  the  matrix  equation  (24)  of  Zhu  et  al.  (1992).  With  a  non-FTE-corrected 
C  thus  computed,  the  cooling  rate  (1)  follows  by  specifying  the  wide-band  normalized  source 
function  column  vector  0  =  0(T)  at  the  band  center  frequency  of  vq  ~675  cm'^ 


OiTk)  = 


B(vo,  Tk) 
B(vo,  200  K) 
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127.38 

exp(^)-l’ 


(5) 


where  B{vq,  T^)  is  the  Planck  function.  The  source  function  is  normalized  in  (5)  by  the  Planck 
function  value  at  200  K  so  that  C,  scaled  by  the  same  factor,  yields  approximate  cooling-rate  con¬ 
tributions  in  its  matrix  elements  [this  normalization  factor  cancels  out  of  (1)].  The  normalization 
temperature  of  200  K  in  (1)  is  scaled  down  from  the  250  K  used  for  Earth  by  Zhu  (1993)  to  reflect 
generally  lower  mean  temperatures  on  Mars. 

2.2.  Computation  and  Validation  of  Mean  Rates 

We  defined  representative  ranges  of  atmospheric  variability  by  averaging  vertical  profiles  of 
temperature,  pressure,  and  CO2  and  O  mixing  ratios  at  a  given  latitude  from  the  Mars  Climate 
Database  (MCD:  Lewis  et  ah,  1999).  Solid  curves  in  Figure  1  show  mean  and  maximum/minimum 
profiles  of  various  fields  at  40°N  under  conditions  of  typical  background  dust  loading  and  moderate 
solar  EUV  activity  (Forget  et  ah,  1999).  The  gray  scale  depicts  the  profile  distribution  density  for 
the  kinetic  temperature  T ,  CO2  and  O  mixing  ratios  and  O/CO2  mixing  ratio. 

The  mean  profiles  in  Figure  1  are  used  as  our  reference  state  for  the  cooling-rate  computations 
and  are  interpolated  to  a  height  grid  Zk  (k  =  1 ...  L)  of  resolution  Az  =  0.5  km  at  L  =  401  levels 
from  the  surface  to  200  km.  The  reference  Curtis  matrix  Co  and  the  various  interpolation  matrices 
in  (2)  were  then  calculated  by  the  correlated-k  distribution  method  (Lacis  and  Oinas,  1991;  Zhu, 
1994).  The  correlated-k  coefficients  and  the  temperature  and  pressure  dependences  were  calcu¬ 
lated  by  line-by-line  integration  using  the  1986  HITRAN  database  (Rothman  et  ah,  1987)  for  20 
different  pressure  levels  (from  0  to  1100  hPa)  and  7  reference  temperatures  from  150  K  to  300  K 
(see  Table  1  of  Zhu  et  ah,  1992).  The  entire  CO2  15-//m  band  (525-825  cm“^)  was  divided  into 
12  equally  spaced  subbands  and  30  Gaussian  quadrature  points  were  used  in  each  subband  to  per¬ 
form  the  integration  (Zhu,  1994).  Since  the  vibrational  partition  function  of  the  CO2  15-/im  band 
is  nearly  constant,  the  coefficient  was  simply  extrapolated  to  100  K  to  cover  the  lower  temperature 
range  encountered  in  the  atmosphere  of  Mars  (see  Figure  lb).  The  resulting  reference  matrices 
and  vectors  were  saved  as  lookup  tables  for  subsequent  use  in  (2). 

In  addition  to  the  standard  modifications  of  physical  constants,  the  following  rate  changes 
were  made  for  Mars.  For  Earth,  Zhu  (1990)  used  deactivation  rates  of  CO2  y2  V-T  emissions  by 
collisions  with  air  molecules  (N2  plus  O2)  and  O  atoms  from  Dickinson  (1984)  and  references 
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therein.  For  the  predominantly  CO2  atmosphere  of  Mars,  the  former  rate  [eq.  (73a)  in  Zhu  (1990)] 
is  replaced  with  a  new  rate  based  on  collisional  deactivation  by  other  CO2  molecules.  For  the 
y2  =  1  — >  0  V-T  transition  we  adopt  the  analytical  fit  {k^a)  in  Table  3  of  Lopez- Valverde  and  Lopez- 
Puertas  (1994).  The  rate  due  to  collisional  deactivation  by  monatomic  oxygen  (O^P),  ko-co2^ 
which  has  a  large  impact  on  high-altitude  cooling  rates,  was  originally  2  x  10"^^  cm  s“^  in  (73b) 
of  Zhu  (1990).  It  is  reset  here  to  3  x  10“^^  cm  s“',  based  on  recent  modeling  and  laboratory  work 
(e.g..  Rougher  et  ah,  1994;  Khvorostovskaya  et  ah,  2002;  Castle  et  ah,  2006),  and  is  identical  to  the 
value  for  Mars  used  by  Lopez- Valverde  and  Lopez-Puertas  (1994)  (their  kja)  and  recommended  by 
Lopez- Valverde  et  al.  (1998).  Lacking  measurements  of  the  rates  of  y2  =  2  — >  1  V-T  transitions, 
Lopez- Valverde  and  Lopez-Puertas  (1994)  scaled  their  y2  =  1  — >  0  rates  by  2.5  and  2  for  collisions 
with  CO2  and  O^P,  respectively,  which  we  also  adopt,  the  former  rate  replacing  (72b)  of  Zhu 
(1990). 

Figure  2a  plots  the  mean  kinetic  temperature  profile  used  by  Lopez- Valverde  et  al.  (1998)  to 
compute  cooling  rates  using  the  more  detailed  Curtis-matrix  Mars  RT  code  of  Lopez- Valverde 
and  Lopez-Puertas  (1994).  Figure  2b  plots  the  15  /rm  cooling  rates  from  our  interpolated  Curtis 
matrix  calculation  using  their  temperature  profile  and  mean  profiles  of  O  and  CO2  mixing  ratio 
in  Figures  Ic  and  le,  respectively.  Despite  differences  in  mixing  ratio  profiles  and  our  simpler 
RT  model,  the  cooling  rates  in  Figure  2b  are  very  similar  in  overall  shape  and  magnitude  to  those 
of  Lopez- Valverde  et  al.  (1998)  (see  their  Figure  lb).  Figure  3c  shows  variations  in  cooling  rate 
for  the  mean  40°N  profiles  in  Figure  1  as  ko-co2  is  halved  and  doubled.  The  cooling  rates  show  a 
similar  range  of  variability  to  that  found  in  a  similar  sensitivity  test  of  Mars  cooling  rates  performed 
by  Lopez-Puertas  and  Lopez- Valverde  (1995)  (see  their  Figure  6).  Other  tests  (not  shown)  reveal 
similarly  strong  increases  to  those  found  by  Lopez-Puertas  and  Lopez- Valverde  (1995)  when  the 
profile  of  O^P  mixing  ratios  is  increased  to  values  representative  of  solar  maximum  conditions. 

Here  and  throughout  this  paper,  we  express  our  heating  rates  and  damping  rates  in  terms  of 
Earth  days  (day“^).  Corresponding  rates  in  terms  of  a  Martian  day  (soL^)  are  obtained  by  multi¬ 
plying  these  rates  by  1.027491. 
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3.  Vertical  Scale-Dependent  IR  Radiative  Damping  Rates 


Our  approach  to  quantifying  scale-dependent  IR  radiative  damping  of  temperature  perturba¬ 
tions  closely  follows  that  of  Zhu  and  Strobel  (1991)  and  Zhu  (1993),  who  provide  more  details 
and  justification. 

If  we  add  small  temperature  perturbations  dT  to  a  background  temperature  vector  Tb,  lineariza¬ 
tion  of  (1)  yields 

dST 
dt 

where  terms  computed  at  the  unperturbed  temperature  Tb  have  the  subscript  “b.”  Eq.  (6)  can  be 
reexpressed,  with  the  aid  of  matrix  manipulations  (Zhu  and  Strobel,  1991),  as 


dT 


6T  +  dC0b, 


(6) 


d6T 

dt 


DdT. 


(7) 


Calculations  by  Zhu  and  Strobel  (1991)  show  that,  for  small  dT,  contributions  from  the  first  term 
in  (6)  are  an  order  of  magnitude  larger  than  from  the  second.  Thus,  to  a  good  approximation 


D  = 


dT 


(8) 


The  linear  system  (7)  and  (8)  describes  nonlocal  damping  of  temperature  perturbations,  and 
can  be  solved  using  either  the  eigenvalue  or  the  Newtonian  cooling  approach  (Goody  and  Yung, 
1989).  Both  solutions  were  computed  and  compared  by  Zhu  and  Strobel  (1991)  and  Zhu  (1993). 
They  found  that,  while  each  gave  equivalent  results  in  certain  limits,  the  two  solutions  generally 
differed.  They  found  that  the  most  accurate  solution  overall  was  a  Newtonian  cooling  rate  that 
varied  both  with  height  and  the  vertical  scale  of  the  temperature  perturbation,  derived  explicitly  as 
(Sasamori  and  London,  1966;  Zhu  and  Strobel,  1991) 


:\z,m)  =  -  f 
%J  — 0> 


K(z',  z)  cos[m(z'  -  z)]dz'. 


(9) 


where  K{z' ,z)  is  the  kernel  function  (Spiegel,  1957),  m  =  In/A^  is  the  vertical  wavenumber  of 
(assumed)  sinusoidal  eigenmodes,  and  Tr{z,  ni)  is  the  e-folding  damping  time  scale  of  a  sinusoidal 
temperature  perturbation  of  wavenumber  m  due  to  infrared  radiative  cooling.  Since  (7)  and  (8)  can 
be  expressed  equivalently  as  a  discretized  form  of  the  integral  relation  (Zhu,  1993) 

— - — =  K(z  ,z)ST(z)dz  ,  (10) 

dt  J_oo 
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then  solutions  to  (9)  follow  by  evaluating  (7)  and  (8)  at  Zi  using  a  perturbation  temperature  profile 
6T(z')  =  t  cos  m(z'  -  Zi),  whereupon  T~^(zi,m)  =  -t~^d6T(zi)/dt.  The  choice  for  f  is  arbitrary 
since  it  has  no  impact  on  the  final  m)  estimate  due  to  the  linear  nature  of  the  calculation. 

We  computed  rates  using  the  mean  temperature  and  constituent  profiles  in  Figure  1  as  our 
background,  for  =  In/m  values  in  the  range  1-500  km.  Figure  4  profiles  the  resulting  t~^(z,  m) 
estimates.  We  see  strong  vertical  scale  dependence  in  the  rates  at  z  <  100  km,  with  damping 
rates  peaking  at  ~30-40  days“^  (i.e.,  e-folding  times  of  <1  hour)  at  the  shortest  vertical  scales 
at  z  ~60-80  km.  The  rates  become  largely  scale-independent  at  z  >120  km  due  to  breakdown  of 
LTE  that  yields  an  optically  thin,  transparent  damping  limit  that  is  well  approximated  by  a  constant 
cooling-to-space  term,  due  to  negligible  layer  exchange. 

Figure  5  profiles  these  radiative  damping  rates  for  different  choices  of  A^  =  In/m.  This  plot 
resembles  the  presentation  in  Figure  4  of  Zhu  (1993)  of  corresponding  CO2  and  O3  damping  rates 
for  the  Earth’s  atmosphere.  Comparisons  between  the  two  plots  show  that  scale-dependent  thermal 
damping  rates  due  to  15  yum  CO2  cooling  in  the  Martian  atmosphere  are,  overall,  about  20  times 
faster  than  on  Earth.  This  finding  is  broadly  consistent  with  an  earlier  Mars-Earth  comparison  of 
radiative  damping  rates  by  Imamura  and  Ogawa  (1995).  Eigure  6  compares  our  scale-dependent 
rates  with  those  of  Imamura  and  Ogawa  (1995),  which  were  derived  using  the  analytical  approx¬ 
imations  of  Eels  (1984).  Our  values  are  generally  larger  than  theirs.  A  similar  offset  between  the 
two  rate  calculations  on  Earth  was  noted  by  Zhu  (1993),  and  ascribed  to  a  more  accurate  treatment 
of  non-ETE  effects  in  the  current  approach  relative  to  that  of  Eels  (1982,  1984).  Given  other  large 
differences  between  the  two  calculations  (e.g.,  different  mean  temperatures,  mixing  ratios  and  rate 
coefficients;  see  below),  there  is  reasonable  overall  agreement  in  Eigure  6  between  the  rates  from 
each  calculation  at  all  altitudes.  Note  in  particular  the  transition  from  strong  scale  dependence  at 
z  =  60-80  km  to  an  essentially  scale-independent  rate  at  z  >  120  km. 

Our  rates  were  computed  using  the  mean  background  temperature  Ti,(z)  in  Eigure  lb,  repro¬ 
duced  as  the  thick  curve  in  Eigure  7a,  with  the  damping  rate  profile  it  yields  at  A^  =5  km  also 
plotted  in  Eigure  7b.  Eollowing  Eels  (1982),  we  can  scale  those  reference  rates  to  yield  damping 
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rates  at  some  other  background  temperature  T  as 


de/dT 


where,  from  (5), 

de  I  exp  (22^)  1 

—  =  1.23686  X  10'  K  — - /  .  (12) 

ir4exp(22p^)-lfj 

Thin  curves  in  Figure  7a  plot  temperature  profiles  offset  from  Ti,{z)  by  +5  K,  +10  K  and  ±20  K. 
The  corresponding  temperature-scaled  rates  m;  T)  at  Tj,  =  5  km  from  (11)  and  (12)  are  plot¬ 
ted  with  thin  curves  in  Figure  7b.  The  rates  show  a  strong  dependence  on  background  temperature 
through  (12).  The  impact  is  particularly  large  at  z  ~70-120  km  due  to  the  low  background  tem¬ 
peratures  at  these  altitudes  which  yields  a  larger  relative  change  to  (12),  such  that  a  uniform  10  K 
warmer  atmosphere  at  all  altitudes  causes  a  doubling  of  the  radiative  damping  rate  at  these  alti¬ 
tudes. 

Appendix  A  describes  least-squares  fits  of  our  Th)  rates  that  can  be  used  with  (11) 

and  (12)  as  a  simple  computationally-efRcient  parameterization  of  scale-dependent  IR  radiative 
damping  rates. 


4.  Gravity  Wave  Damping 

Following  Marks  and  Eckermann  (1995),  we  consider  gravity  waves  governed  by  the  anelastic 
dispersion  relation 

m  = - - — - - a  .  (13) 

-  r 

Wave  parameters  are  as  follows:  oj  is  ground-based  frequency,  co^  =  co  -  Uk  -  Vlis  intrinsic  fre¬ 
quency,  and  (k,  I,  m)  is  the  wavenumber  vector.  Background  atmospheric  variables  are  horizontal 
wind  (t/,  V),  inertial  frequency  /,  buoyancy  frequency  N,  and  =  1  where  Hp  is  the  density 
scale-height. 

Gravity  wave  dissipation  affects  the  larger-scale  circulation  through  deposition  of  the  wave 


momentum  flux 


F  =  pfo  {u'w' ,  v'w'^  =  (k,  1)  CgzA. 
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Here  (u',  V ,  w')  is  the  gravity  wave-induced  vector  velocity  anomaly,  pb  is  background  density,  Cg^ 
is  vertical  group  velocity,  A  =  E  Ico'^  is  wave  action  density,  E  is  total  (kinetic  plus  potential)  wave 
energy,  and  overbars  denote  horizontal  averages. 

In  the  absence  of  dissipation,  wave  action  density  A  is  conserved  through  the  continuity  relation 

dA 

— +V.(CgA)  =  0,  (15) 

where  Cg  is  the  ground-based  three-dimensional  group  velocity  vector.  In  the  absence  of  horizontal 
gradients  and  local  time  variations,  the  solution  of  (15)  is  a  constant  vertical  flux  of  wave  action, 
Cg^A,  and  thus  constant  momentum  flux  F  through  (14),  since  (k, /)  is  also  constant.  Dissipation 
of  gravity  waves  causes  some  fraction  of  the  wave  momentum  flux  (14)  to  be  deposited  locally, 
which  in  turn  exerts  a  force  -i)F  jdz  on  the  background  flow.  Appendix  B  of  Marks  and  Eckermann 
(1995)  outlines  how  general  mechanical  and  thermal  dissipation  terms  in  the  linearized  Navier 
Stokes  equations  yield  net  wave  action  dissipation  terms  on  the  right  hand  side  of  (15).  Substituting 
standard  expressions  for  viscous  dissipation  and  thermal  Newtonian  damping,  the  net  dissipation 
of  wave  action  density  can  be  expressed  in  the  Newtonian  damping  form  (Marks  and  Eckermann, 
1995) 


dA  2A 

^^+V.(c,A)- 


where  the  damping  rate  of  gravity-wave  amplitudes 


+^-i  l  +  ^  + 


\-fl0J*  I  Pj.-1 


\  +  ^  + 


Twr  +  ^ 


N^-(o+ 

-1 


(16) 


(17) 

(18) 


In  (17),  is  the  scale-dependent  IR  radiative  damping  rate  (11),  =  ri(k^  +  f  +  np-  +  a^) 

is  the  viscous  damping  rate,  77  is  the  kinematic  viscosity  and  Pr  is  the  Prandtl  number.  In  the 
equivalent  form  (18),  and  t~1  are  the  damping  rates  of  gravity-wave  amplitudes  due  to  IR 
radiative  cooling  and  viscous  dissipation,  respectively.  These  in  turn  depend  on  the  gravity-wave 
weighting  coefflcents  for  radiative  damping,  Wr,  and  viscous  damping,  Wy.  Eigure  8  shows  how 
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they  vary  as  intrinsic  frequency  co""  is  varied  from  the  low-frequency  limit  {f^  jco*^  ~  1)  through 
to  a  nominal  high-frequency  limit  of  1  (note  that  gravity  waves  governed  by  (13)  have 

a  high-frequency  limit  oj*  that  is  typically  less  than  N\  see  Marks  and  Eckermann,  1995).  The 
variations  in  Figure  8  are  due  to  the  polarization  characteristics  of  gravity  waves.  For  example, 
the  vanishing  of  Wy  (and  hence  t“J)  as  is  due  to  an  enhancement  in  wave  kinetic  energy 

at  the  expense  of  wave  potential  energy,  which  reduces  wave-induced  temperature  perturbations 
and  thus  the  net  impact  of  thermal  damping  on  wave  action.  Note  too  the  strong  sensitivity  of  the 
viscous  damping  of  gravity  waves  to  the  Prandtl  number. 

Below  the  turbopause,  viscous  damping  has  both  molecular  and  turbulent  contributions.  Since 
the  Prandtl  numbers  of  each  vary,  to  treat  both  forms  of  viscous  damping  it  proves  necessary  to 
further  subdivide  the  viscous  gravity- wave  damping  in  (18)  into  the  sum  of  separate  molec¬ 
ular  and  turbulent  contributions,  given  in  each  case  by  the  second  term  on  the  right-hand  side 
of  (17),  but  in  each  instance  substituting  the  relevant  molecular  or  turbulent  viscosity  (rj^  or  77,, 
respectively)  and  Prandtl  number  (Pr^  or  Pvt,  respectively):  i.e., 

=  W^mTvl  +  (19) 

Since  turbulent  viscosities  rjt  on  Mars  are  poorly  constrained  (see,  e.g.,  Bittner  and  Fricke, 
1987;  Chassefiere  et  ah,  1994;  Izakov,  2007),  we  defer  their  treatment  and  discussion  until  later 
in  the  paper,  and  consider  for  now  just  the  molecular  viscosity  contribution.  We  specify  kinematic 
molecular  viscosities  as 

9.18X  lo-^r^' 

nm  =  - (20) 

Pb 

The  numerator  in  (20)  is  a  least-squares  fit  to  molecular  viscosities  (in  units  of  kg  m“^  s“^)  over 
the  50-500  K  range  from  laboratory  measurements  for  a  Mars-like  mixture  of  gases,  taken  from 
Table  1  of  Catalfamo  et  al.  (2009).  In  the  r^=100-400  K  range,  data  in  Tables  1-2  of  Catalfamo 
et  al.  (2009)  yield  Pr^  ~0.7-0.95,  while  a  representative  Martian  Cp  =  805  J  kg“^  K“^  yields 
Pr^  ~0.89  when  applied  to  the  same  data.  We  adopt  a  median  value  of  Pr^  =  0.8,  slightly  larger 
than  Pr^  =  0.71  for  the  Earth’s  atmosphere,  consistent  with  theory  and  laboratory  measurements 
that  a  mostly  triatomic  molecular  atmosphere  of  CO2  should  have  a  slightly  higher  Prandtl  number 
than  the  mostly  diatomic  atmosphere  of  Earth. 
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The  inverse  density  dependence  in  (20)  yields  an  exponential  increase  in  rj^  with  height,  so 
that  local  molecular  viscous  damping  will  eventually  exceed  the  local  IR  radiative  damping  of 
waves  at  some  altitude  Zy.  Imamura  and  Ogawa  (1995)  estimated  Zy  to  be  ~120  km  on  Mars. 
We  used  the  gravity-wave  equations  above  to  calculate  Zy  based  on  mean  temperature  and  density 
profiles  in  Figure  1  that  specify  both  our  rates  from  Figure  4  and  our  values  from  (20).  The 
profiles  of  N  and  Hp  based  on  this  mean  temperature  profile  are  plotted  in  Figure  9.  To  keep  these 
initial  calculations  simple,  here  we  adopt  constant  values,  based  on  approximate  profile  means  in 
Figure  9,  of  N  =0.01  rad  s“^  and  Hp  =9.5  km,  so  that  m  in  (13)  remains  constant  (later  calculations 
will  remove  this  assumption).  We  set  /  =10“"^  s“',  corresponding  to  a  latitude  of  45°N.  Figure  10 
plots  the  altitude  at  which  first  exceeds  t~].  as  a  function  of  vertical  wavenumber  m  and  total 
horizontal  wavenumber  Ktot  =  {k^  +  We  find  a  strong  variation  of  Zv  with  m,  varying  from 

~80  km  at  Tj;  ~1  km  to  ~140  km  at  ~100  km. 

To  study  in  greater  depth  the  impact  of  radiative  damping  relative  to  other  gravity-wave  damp¬ 
ing  mechanisms  on  Mars,  we  developed  a  simple  one-dimensional  anelastic  gravity- wave  model 
based  on  the  dispersion  relation  (13)  and  it’s  associated  linear  gravity-wave  equations,  as  described 
in  Appendix  B.  The  model  simulates  gravity- wave  propagation  and  amplitude  evolution  through 
a  vertical  atmospheric  column  extending  from  the  surface  to  200  km  using  the  same  height  grid 
Zk  as  for  the  radiation  calculations.  Peak  vertical  displacement  amplitudes  ^(z)  and  momentum 
fluxes  F(z)  =  |F(z)|  are  computed  subject  to  radiative  and  molecular  viscous  damping,  and  am¬ 
plitude  saturation  due  to  wave  breaking.  The  model  is  a  slightly  more  complicated  version  of 
the  one-dimensional  formulations  used  in  gravity-wave  drag  parameterizations  (Joshi  et  ah,  1996; 
Collins  et  ah,  1997;  Forget  et  ah,  1999).  Although  background  horizontal  winds  U(z)  and  y(z) 
are  included  in  the  model,  for  simplicitly  we  set  [/  =  V  =  0  so  that  large  wavenumber  refraction 
effects  are  eliminated,  allowing  us  to  more  effectively  study  dissipation  as  a  function  of  vertical 
wavenumber  at  a  range  of  altitudes.  Since  IR  radiative  damping  rates  are  acutely  temperature 
sensitive  (see  Figure  7),  the  simulated  gravity  waves  propagate  through  the  realistic  background 
temperature  Ttiz)  shown  in  Figures  lb  and  7a.  For  consistency,  our  model  calculations  include  the 
resulting  height  variations  in  N  and  Hp  shown  in  Figure  9,  which  cause  m  values  to  vary  slightly 
with  height  through  the  dispersion  relation  (13). 
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For  each  wave,  we  must  choose  a  launch  (source)  altitude,  and  an  initial  peak  vertical  dis¬ 
placement  amplitude,  i(zs).  We  choose  =  0  on  the  assumption  that  major  gravity-wave  sources 
are  located  in  the  lower  atmosphere,  such  as  waves  forced  by  flow  over  mesoscale  terrain  (Joshi 
et  ah,  1996).  There  are  few  observations  to  reliably  guide  our  choice  for  ^(zj.  On  Earth,  mean 
gravity-wave  amplitudes  are  controlled  by  vertical  wavenumber  spectra  with  highly  reproducible 
shapes  that  peak  in  variance  at  a  characteristic  vertical  wavenumber,  m*,  which  is  typically  in  the 
range  2;r(l-3  km)"'  in  the  lower  atmosphere  (Allen  and  Vincent,  1995).  While  there  have  been  no 
measurements  of  vertical  wavenumber  spectra  on  Mars  to  date,  the  observed  horizontal  wavenum¬ 
ber  spectra  of  Imamura  et  al.  (2007)  bear  a  strong  resemblance  to  those  measured  on  Earth  in  the 
mesoscale  range,  which  in  turn  suggests  thoeretically  that  vertical  wavenumber  spectra  should  also 
have  somewhat  similar  shapes  (see,  e.g.,  Bacmeister  et  al.,  1996).  The  total  vertical  displacement 
variance,  obtained  by  integrating  canonical  forms  for  gravity-wave  vertical  wavenumber  spectra 
observed  on  Earth,  is  ~  O.lmJ  (Eritts  and  VanZandt,  1993),  which,  on  choosing  a  representative 
m*  =  2;r(2  km)"'  for  the  lower  atmosphere,  yields  an  r.m.s.  wave-induced  vertical  displacement  of 
100  m.  Given  conflicting  evidence  as  to  whether  gravity-wave  amplitudes  are  weaker  or  stronger 
on  Mars  relative  to  Earth  (see  section  5),  we  adopt  the  default  assumption  of  similar  wave  ampli¬ 
tudes  and  thus  set  ^(zj  =100  m  for  our  simulations. 

Eigure  1 1  plots  ^(z)  (left  panel)  and  F(z)  (right  panel)  for  a  gravity  wave  of  T/,  =  InjKtot  = 
100  km  and  ground-based  phase  speed  c  ~3.6  m  s“'  which  yields  a  source-level  vertical  wave¬ 
length  Tj,  =  27r/|m(Zi)|  =  2.5  km.  The  thick  gray  curves  in  Eigure  11  show  gravity-wave  solutions 
in  which  no  dissipation  has  been  applied,  whereupon  F{z)  remains  constant  with  height  and  ^(z) 
grows  with  height  roughly  as  oc  exp  f  dzl2Hp.  The  thick  gold  curves  are  solutions  in  which 
only  wave  breaking/saturation  has  been  applied,  which  is  the  standard  parameterization  approach 
(Joshi  et  al.,  1996;  Collins  et  al.,  1997;  Eorget  et  al.,  1999).  The  black  dashed  (blue  dotted)  curve 
is  the  solution  in  which  only  IR  radiative  (molecular  viscous)  damping  has  been  applied,  with 
the  solid  aqua  curve  showing  solutions  with  both  IR  radiative  and  molecular  viscous  damping 
applied.  The  red  curve  shows  the  solution  in  which  all  three  damping  processes  (IR  radiative 
damping,  molecular  viscous  damping,  and  saturation  due  to  wave  breaking)  were  applied.  The 
results  in  Eigure  1 1  show  that  the  amplitude  evolution  of  this  wave  is  dominated  by  IR  radiative 
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damping.  The  IR  damping  of  this  wave  is  so  severe  that  ^(z)  decreases  with  altitude  in  Figure  1  la 
and  never  attains  a  threshold  amplitude  for  wave  breaking. 

Figure  12  shows  corresponding  solutions  for  a  wave  of  ground-based  phase  speed  c  ~6.6  m  s”\ 
which  yields  a  source-level  vertical  wavelength  =  2nl\m{zs)\  =5  km.  While  wave  amplitudes 
grow  slowly  with  altitude  in  this  example  up  to  ~40  km  altitude,  they  are  strongly  dissipated  above 
40  km  and  again  do  not  attain  breaking  amplitudes  at  any  altitude.  Dissipation  of  wave  amplitude 
and  momentum  flux  is  again  dominated  by  the  IR  radiative  damping  term. 

Figure  13  plots  solutions  for  a  wave  of  c  ~12.9  m  s“\  which  yields  a  source-level  verti¬ 
cal  wavelength  A^  =  2nl\mizs)\  =10  km.  This  wave  grows  in  amplitude  with  height  and  attains  a 
breaking  amplitude  at  ~60  km,  whereupon  it  continues  breaking  as  a  saturated  wave  up  to  ~  100  km 
before  being  rapidly  dissipated  by  molecular  viscosity  thereafter.  Clearly,  all  three  damping  mech¬ 
anisms  (IR  radiative  cooling,  molecular  viscosity,  and  wave  breaking)  play  a  role  in  progressively 
depositing  this  wave’s  momentum  flux  throughout  the  atmospheric  column.  To  assess  the  relative 
contributions,  the  red  curve  in  Figure  13c  plots  the  cumulative  percentage  of  the  wave’s  original 
momentum  flux  that  has  been  deposited  into  the  background  flow  for  the  simulation  in  which  all 
three  damping  processes  were  activated,  showing  that  essentially  100%  of  the  wave’s  momentum 
flux  has  been  dissipated  by  the  time  the  wave  has  propagated  to  z  ~  100  km.  The  remaining 
curves  in  Figure  13c  show  the  contributions  to  this  total  from  saturation  (gold  curve),  IR  radiative 
damping  (black  dashed  curve),  molecular  viscosity  (blue  dotted  curve)  and  radiative  plus  molec¬ 
ular  damping  (aqua  curve).  These  results  show  that  ~90%  of  this  wave’s  momentum  flux  was 
dissipated  by  IR  radiative  damping,  ~I0%  by  wave  breaking  (saturation),  and  <1%  by  molecular 
viscosity. 

Only  on  progressing  to  Figure  14,  which  plots  solutions  for  a  wave  of  c  ~19.0  m  s“^  (A^  = 
2n/\m(zs)\  =15  km),  do  we  see  wave  breaking/saturation  become  competitive  with  IR  radiative 
damping  in  dissipating  wave  momentum  flux.  In  this  example.  Figure  14c  reveals  that  wave  break¬ 
ing  accounts  for  slightly  more  than  50%  of  the  overall  flux  deposition  through  the  column,  and 
IR  radiative  damping  accounts  for  slightly  less  than  50%.  Again,  molecular  viscous  damping  is 
essentially  negligible,  acting  only  to  damp  out  the  tiny  momentum  flux  residual  that  remains  at 
high  altitudes. 
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To  provide  a  broader  integrated  picture,  we  conducted  a  series  of  simulations  like  those  in 
Figures  11-14  for  a  range  of  waves  with  different  source-level  vertical  wavenumbers  m(zs),  with 
all  physical  damping  processes  activated  in  each  case.  Figure  15a  plots  the  cumulative  percentage 
of  dissipated  momentum  flux  as  a  function  of  z  and  |m|  for  waves  of  Ktot  =  2n{l00  km)“^  The 
contributions  to  the  percentage  totals  in  Figure  15a  from  wave  breaking/saturation,  molecular 
viscous  damping,  and  IR  radiative  damping,  are  plotted  in  Figures  15b,  15c  and  15d,  respectively. 
In  all  cases,  molecular  viscous  damping  accounts  for  much  less  than  1%  of  the  total  flux  deposition. 
At  low  altitudes,  IR  radiative  damping  is  the  dominant  momentum-flux  dissipation  mechanism  for 
all  the  waves,  and  it  remains  dominant  at  higher  altitudes  for  those  gravity  waves  of  /Iz(Zi)  < 
15  km.  At  higher  altitudes,  saturation  eventually  becomes  the  dominant  flux  deposition  process 
for  gravity  waves  of  T^CZi)  ^15  km,  accounting  for  >  90%  of  the  flux  deposition  for  waves  with 
>  40  km. 

These  findings  are  highlighted  in  a  slightly  different  way  in  Figures  15e  and  15f,  which  plot 
the  ratio  of  saturation  and  IR  radiative  damping  percentages,  respectively,  to  the  overall  local 
flux  deposition  percentage  in  Figure  15a.  Contours  exist  only  in  the  top-left  regions  of  these 
panels,  with  the  remaining  “blank”  regions  corresponding  to  a  broad  area  of  m-z  space  where  the 
saturation  contribution  to  the  momentum-flux  dissipation  is  <1%,  and  the  IR  radiative  damping 
contribution  is  >99%.  This  highlights  the  major  role  played  by  IR  radiative  damping  in  controlling 
the  momentum  flux  deposition  of  these  waves. 

Figure  16  plots  “top-of-the-atmosphere”  momentum-flux  deposition  percentages  at  z  =  200  km 
due  to  each  gravity-wave  dissipation  process,  resulting  from  3600  individual  gravity-wave  simu¬ 
lations  spanning  the  indicated  range  of  horizontal  and  vertical  wavenumber  pairings.  Figure  16b 
again  reveals  a  large  area  of  this  two-dimensional  wavenumber  space  where  IR  radiative  damping 
is  the  dominant  flux  deposition  process.  Molecular  viscous  damping  is  again  very  small  for  every 
wavenumber  pairing  considered.  Wave  breaking  is  significant  only  at  longer  vertical  wavelengths 
and  shorter  horizontal  wavelengths,  which  is  consistent  with  the  faster  vertical  group  velocities 
at  these  wave  scales  that  yield  less  time-integrated  IR  radiative  damping  along  the  ray  path,  in 
addition  to  the  smaller  radiative  damping  rates  that  occur  at  long  vertical  wavelengths  in  Figure  4. 
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5.  Discussion  and  Conclusions 


While  faster  radiative  relaxation  rates  on  Mars  relative  to  Earth  were  inferred  more  than  40 
years  ago  by  Goody  and  Belton  (1967),  the  damping  of  small  vertical-scale  temperature  structure 
on  Mars  by  15  //m  CO2  cooling  has  remained  largely  unquantified.  In  particular,  it  has  remained 
unclear  whether  longwave  radiative  cooling  yields  significant  thermal  damping  of  gravity-wave 
momentum  fluxes.  The  parameterizations  of  gravity-wave  drag  currently  included  within  Mars 
GCMs  do  not  include  it,  and  thus  tacitly  assume  that  it  is  a  secondary,  negligible  effect. 

Using  a  Curtis  matrix  model  of  15  /rm  radiative  cooling  specifically  adapted  for  Mars,  we 
computed  scale-dependent  IR  radiative  damping  rates  ni)  from  z  =0-200  km  over  a  vertical 
wavelength  band  from  1-500  km,  for  a  representative  mean  Martian  atmosphere  at  40°N.  We  then 
used  those  rates  to  study  the  impact  of  IR  radiative  damping  on  the  dissipation  of  gravity  wave 
momentum  fluxes,  and  compared  it  to  those  of  other  gravity-wave  dissipation  processes,  such  as 
wave  breaking  (saturation)  and  molecular  viscous  dissipation.  Using  a  linear  gravity  wave  model 
governed  by  simple  anelastic  wave  physics,  our  results  suggest  that  IR  radiative  damping  is  the 
dominant  process  by  which  gravity  wave  momentum  flux  is  dissipated  on  Mars  at  altitudes  below 
~50  km,  and  is  also  the  dominant  flux-dissipating  process  at  all  altitudes  for  those  waves  whose 
vertical  wavelengths  are  <10-15  km  (see  Figures  15  and  16). 

The  precise  demarcations  with  altitude  and  wavenumber  summarized  in  Figures  15  and  16 
depend  on  a  number  of  specific  assumptions  adopted  in  our  modeling.  For  example,  source-level 
gravity-wave  amplitudes  t{Zs)  could  be  much  larger  than  our  current  choice  of  100  m,  which  would 
cause  waves  to  break  at  lower  altitudes  and  thus  increase  the  role  of  wave  breaking.  Different 
background  thermal  conditions  will  also  reduce  (or  increase)  radiative  damping  rates.  We  ran  a 
series  of  simulations  in  which  these  and  other  choices  were  varied.  While  the  locations  of  the 
contours  in  Figures  15  and  16  can  vary  in  response,  the  overall  morphology  does  not  change 
markedly,  and  hence  the  major  conclusions  that  we  drew  from  the  large-scale  features  evident  in 
these  results  also  do  not  change  in  any  substantial  way:  Specifically, 

•  IR  radiative  damping  is  the  dominant  momentum-flux  dissipation  mechanism  for  all  gravity 

waves  in  the  lower  atmosphere  and  lower  regions  of  the  middle  atmosphere; 
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•  for  gravity  waves  whose  souree-level  vertieal  wavelengths  are  <10-15  km,  IR  radiative 
damping  is  the  dominant  momentum-flux  dissipation  meehanism  at  all  altitudes; 

•  moleeular  viseosity  has  negligible  impaet  on  the  total  flux  deposition  of  gravity  waves;  and, 

•  wave  breaking  is  the  dominant  flux  dissipation  meehanism  for  “fast”  gravity  waves  of  short 
horizontal  wavelength  and/or  long  vertieal  wavelengths,  whieh  eseape  major  IR  damping  at 
lower  altitudes  and  deposit  the  majority  of  their  momentum  flux  through  wavebreaking  and 
saturation  at  higher  altitudes. 

In  short,  our  results  strongly  suggest  that  seale-dependent  IR  radiative  damping  is  a  first-order, 
and  frequently  dominant,  meehanism  by  whieh  gravity- wave  amplitudes  are  dissipated  on  Mars, 
and  henee  exerts  an  important  influenee  on  gravity-wave  momentum  flux  deposition  and  the  wave- 
driven  eireulation  of  Mars. 

Parameterizations  of  subgrid-seale  orographie  gravity-wave  drag  eurrently  embedded  within 
Mars  GCMs  are  adapted  from  similar  parameterizations  in  Earth  GCMs,  whieh  only  parameterize 
flux  deposition  due  to  wave  breaking  and  saturation  (Joshi  et  ah,  1996;  Collins  et  ah,  1997;  Forget 
et  ah,  1999).  As  Read  and  Lewis  (2004)  diseuss  in  their  review  of  gravity-wave  drag  etfeets  on 
the  Martian  elimate  (ehapter  4.4.4),  when  these  parameterizations  were  initially  adapted  for  and 
embedded  within  Mars  GCMs,  the  parameterized  souree-level  wave  momentum  fluxes  had  to  be 
redueed  by  1-2  orders  of  magnitude  to  produee  a  realistie  elimate  (e.g.,  Collins  et  ah,  1997).  Taken 
at  faee  value,  this  would  imply  mueh  weaker  orographie  foreing  of  gravity-wave  momentum  flux 
on  Mars  relative  to  Earth.  Yet  surfaee  flow  patterns,  and  partieularly  orographie  relief,  would 
(if  anything)  seem  to  suggest  the  opposite,  and  the  limited  direet  measurements  of  high-altitude 
gravity  waves  on  Mars  to  date  seem  to  support  this  eountervailing  viewpoint  (e.g.,  Fritts  et  ah, 
2006). 

The  ealeulations  presented  in  this  paper  lead  us  to  propose  an  alternative  explanation  of  this 
result.  Our  results  reveal  appreeiable  IR  radiative  damping  of  the  momentum  flux  of  all  gravity 
wave  at  altitudes  z  ~0-50  km.  This  radiative  damping  prevents  short  vertieal  wavelength  gravity 
waves  from  attaining  breaking  amplitudes  (e.g..  Figures  11  and  12),  and  eauses  waves  of  inter¬ 
mediate  and  longer  vertieal  wavelength  to  break  at  higher  altitudes  (e.g..  Figure  13),  where  the 
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wave  possesses  less  remaining  momentum  flux  to  deposit  via  wave  breaking  and  thus  produces  a 
weaker  local  mean-flow  acceleration.  For  parameterizations  that  lack  a  model  of  scale-dependent 
IR  radiative  damping,  the  only  way  to  tune  the  parameterization  to  compensate  for  this  missing 
physics  is  to  articificially  reduce  wave  amplitudes  at  the  source,  so  that  waves  break  higher  and 
deposit  less  momentum  flux  via  breaking. 

More  work  is  needed  to  ascertain  the  exact  ways  in  which  IR  radiative  damping  and  wave 
breaking  combine  to  control  the  deposition  of  gravity-wave  momentum  flux  throughout  the  Mar¬ 
tian  atmosphere  and  in  turn  drive  a  large-scale  wave-driven  circulation.  As  on  Earth,  a  major 
controlling  factor  will  be  the  form  of  the  background  Martian  winds  Uiz)  and  Viz).  In  this  initial 
work,  we  deliberately  omitted  the  complicating  effects  of  mean  wind  shear  and  vertical  wavenum¬ 
ber  refraction  by  setting  U  =  V  =  0,  to  focus  more  easily  on  the  relative  impacts  of  IR  damping 
relative  to  other  processes  as  a  function  of  height  and  wave  scale,  using  waves  of  roughly  con¬ 
stant  frequency  and  wavenumber.  Traditionally,  westerly  zonal  winds  Uiz)  and  westerly  shear, 
for  example,  are  viewed  as  progressively  removing  gravity  waves  with  westerly  phase  speeds  c  at 
critical  levels  as  — >  /^.  Our  results  strongly  suggest  that  such  waves  will  be  heavily  dissipated 

by  IR  radiative  damping  as  their  vertical  wavelengths  contract,  and  will  be  “radiatively  quenched” 
well  before  reaching  either  wave-breaking  amplitudes  or  being  absorbed  at  a  critical  level.  Those 
easterly  phase-speed  waves  propagating  against  the  westerly  flow  will  refract  to  longer  vertical 
wavelengths,  attain  fast  vertical  group  speeds,  and  thus  may  refract  into  regions  of  wavenumber 
space  where  upper-level  wave  breaking  may  be  dominant  (see  Figures  15e  and  16a).  The  situation 
becomes  complicated  since  wave  breaking  conditions  also  vary  strongly  in  response  to  background 
winds  and  shear.  Generalization  to  three  dimensions  with  variable  wind  and  wave  vector  azimuths 
further  complicates  the  picture,  as  does  the  nature  and  distribution  of  wave  sources,  making  es¬ 
timates  of  relative  impacts  complicated  and  case-specific.  While  these  arguments  indicate  that 
more  detailed  model  studies  are  needed  to  fully  quantify  such  effects,  we  believe  that  the  results 
presented  here  already  present  a  sufficiently  strong  and  compelling  case  to  consider  IR  radiative 
damping  an  indispensable  component  of  any  parameterization  that  seeks  to  model  gravity-wave 
momentum  flux  deposition  accurately,  as  is  required  in  GCMs.  As  described  in  Appendices  A 
and  B,  the  first  order  effects  of  IR  radiative  damping  could  be  incorporated  into  existing  gravity- 
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wave  models  and  parameterizations  with  little  increase  in  overall  computational  complexity  and 
overhead. 

To  date,  strong  similarities  have  been  noted  between  gravity-wave  dynamics  on  Earth  and 
Mars  (e.g.,  Pirraglia,  1976;  Pickersgill  and  Hunt,  1979;  Tobie  et  al.,  2003;  Imamura  et  al.,  2007; 
Parrish  et  al.,  2009).  This  is  reflected  in  the  GCM  parameterizations  of  gravity-wave  drag  for  both 
planetary  atmospheres,  which  at  present  possess  almost  identical  governing  physics  (e.g..  Forget 
et  al.,  1999;  Medvedev  and  Hartogh,  2007).  Our  results  point  to  a  major  physical  difference  in 
gravity-wave  dissipation  on  each  planet.  On  Earth,  IR  radiative  damping  is  much  weaker  and  wave 
breaking  is  the  dominant  flux  dissipation  process  (Marks  and  Eckermann,  1995),  and  so  gravity- 
wave  dissipation  is  characterized  by  instability,  overturning,  turbulence  generation  and  turbulent 
mixing.  By  contrast,  we  find  IR  radiative  damping  of  gravity  waves  to  be  a  major  (often  dominant) 
gravity-wave  dissipation  process  on  Mars.  The  IR  radiative  damping  of  gravity-wave  momentum 
and  energy  fluxes,  while  producing  mean- flow  accelerations  and  frictional  heating,  respectively, 
does  not  (in  and  of  itself)  generate  wave-field  instabilities  and  turbulence.  This  in  turn  impacts 
the  potential  role  of  a  fourth  wave-damping  mechanism,  raised  earlier  in  this  paper  [see  eq.  (19)] 
but  thereafter  ignored  in  our  modeling:  dissipation  due  to  background  turbulent  diffusion.  In  the 
Earth’s  middle  atmosphere,  most  turbulence  is  generated  by  gravity-wave  breaking.  If,  as  our 
results  suggest,  IR  radiative  damping  reduces  the  amount  of  gravity-wave  flux  on  Mars  available 
for  wave  breaking  (and  hence  turbulence  generation),  it  would  (all  other  things  being  equal)  reduce 
gravity-wave-induced  turbulence  levels  on  Mars  relative  to  Earth,  as  well  as  the  overall  impact 
of  wave  damping  due  to  background  turbulent  diffusion  and  the  turbulent  vertical  diffusion  and 
mixing  of  trace  constituents  in  the  Martian  middle  atmosphere. 
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Figure  1:  Profiles  of  (a)  pressure,  (b)  temperature,  (c)  CO2  volume  mixing  ratio,  (d)  density,  (e)  O^P  volume  mixing 
ratio,  and  (f)  O^P/C02  volume  mixing  ratio.  Thick  curves  are  the  means  of  individual  profiles  at  40°N  from  the 
Mars  Climate  Database  (MCD:  Lewis  et  al.,  1999)  for  average  levels  of  dust-loading  and  solar  activity.  Thin  curves 
are  minimum  and  maximum  profiles,  and  gray  shading  in  panels  b,  c,  e  and  f  depicts  the  distribution  density  of  the 
individual  MCD  profiles. 
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Figure  2:  (a)  Mean  kinetic  temperature  profile  used  by  Lopez- Valverde  et  al.  (1998)  (their  Figure  la);  (b)  resulting 
15  pm  cooling  rates  using  our  Curtis-matrix  algorithm  and  the  mean  O  and  CO2  mixing  ratio  profiles  from  Figures  Ic 
and  le.  Rates  are  quantified  here  and  elsewhere  in  terms  of  Earth  days:  rates  in  K  sol  '  are  obtained  by  multiplying 
the  values  in  (b)  by  1.027491. 
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Figure  3:  Profiles  of  mean  (a)  kinetic  temperature  and  (b)  O/CO2  mixing  ratio  from  Figures  lb  and  If,  respectively, 
which  yield  in  (c)  a  mean  IR  cooling  rate,  plotted  with  a  thick  curve,  for  our  prespecified  O-CO2  deactivation  rate, 
ko-co2^  of  3  X  10“^^  cm  s  '.  Thin  curves  in  (c)  show  corresponding  cooling  rates  when  this  ko-co2  value  is  doubled 
and  halved. 
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Figure  5:  Scale-dependent  IR  radiative  damping  rates,  t^^(z,  m),  plotted  as  a  function  of  height  z  for  various 
2n/m  values  shown  in  top  right  of  each  panel. 
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Figure  6;  A  comparison  of  IR  radiative  damping  rates,  hi),  at  z  =60,  80,  100,  and  120  km  in  the  present  work 
(solid  curves)  and  corresponding  estimates  of  Imamura  and  Ogawa  (1995)  (dotted  curves). 
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Figure  7:  Variation  of  IR  radiative  damping  rate  m)  at  -  Injm  -  5  km  as  the  reference  temperature  prohle 
in  (a)  (thick  curve)  is  modihed  at  all  heights  by  -20  K,  -10  K,  -5  K,  +5  K,  +10  K  and  +20  K  (thin  curves). 
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Figure  8;  Weighting  factors  Wr  and  Wy  in  (17)  and  (18)  that  scale  the  raw  radiative  and  viscous  damping  rates  '  and 
T“'  to  yield  the  net  rate  of  damping  of  gravity-wave  amplitudes. 
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Figure  9:  Profiles  of  mean  (a)  buoyancy  frequency  N  and  (b)  density  scale  height  Hp,  based  on  the  mean  background 
kinetic  temperature  profile  Thiz)  (thick  solid  curves  in  Figures  lb  and  7a). 
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Figure  1 1 :  Profiles  of  (a)  peak  vertical  displacement  amplitude  ((z)  and  (b)  vertical  flux  of  horizontal  momentum  F{z) 
for  a  gravity  wave  of  c  =  3.6  m  s“',  /l,fc)  =  2jil\m{zs)\  -  2.5  km,  Ak  -  2nlK,ot  -  100  km,  and  l{Zs)  -  100  m.  The 
curves  show  solutions  for  different  types  of  wave  damping;  no  damping  (gray  solid  curves),  wave  breaking/saturation 
only  (gold  solid  curves),  molecular  viscous  damping  only  (blue  dotted  curves),  IR  radiative  damping  only  (black 
dashed  curves),  molecular  viscous  and  IR  radiative  damping  (aqua  solid  curves),  and  all  three  damping  processes 
(wave  breaking/saturation,  molecular  viscous  damping  and  IR  radiative  damping:  red  solid  curves). 


Figure  12: 


As  for  Figure  1 1  but  for  a  gravity  wave  of  c  =  6.6  m  s  '  and  /Ij(zj)  =  2nl\m{zs)\  -  5  km. 
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Figure  13:  Panels  (a)  and  (b)  follow  the  same  presentation  as  Figure  1 1  but  for  a  gravity  wave  of  c  =  12.9  m  s“'  and 
=  27r/|OT(zs)|  =  10  km.  The  red  curve  in  (c)  plots  the  cumulative  deposition  of  momentum  flux  F(z),  expressed 
as  a  percentage  of  the  original  source-level  flux,  for  the  simulation  with  all  three  dissipation  processes  activated  (red 
curves  in  panels  a  and  b).  Remaining  curves  in  (c)  show  the  contributions  to  this  cumulative  flux  deposition  percentage 
from  breaking/saturation  (gold  solid  curve),  IR  radiative  damping  (black  dashed  curve),  molecular  viscous  damping 
(blue  dotted  curve),  and  combined  molecular  viscous  and  IR  radiative  damping  (aqua  solid  curve). 


Figure  14:  As  for  Figure  13  but  for  a  gravity  wave  of  c  =  19.0  m  s  ^  and  A^(zs)  -  27r/|m(zj)|  =  15  km. 
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Figure  15:  Cumulative  percentages  of  dissipated  gravity-wave  momentum  flux  (0%  =  no  dissipation,  100%  =  total 
dissipation)  as  a  function  of  height  z  and  source-level  vertical  wavenumber  \m(zs)\,  for  model-simulated  gravity  waves 
of  K,ot  =  2;7r(100  km)“'  and  ^(Zs)  -  100  m  with  all  damping  mechanisms  activated.  Plots  show  (a)  the  total  percentage 
of  dissipated  flux,  and  then  contributions  to  this  total  from  (b)  breaking/saturation,  (c)  molecular  viscous  damping 
and  (d)  IR  radiative  damping.  Remaining  panels  show  the  percentage  contribution  to  the  local  total  in  (a)  from  (e) 
breaking/saturation  and  (f)  IR  radiative  damping,  the  latter  showing  the  dominant  (i.e.  >99%)  role  of  IR  radiative 
damping  in  gravity-wave  momentum  flux  deposition  for  all  waves  at  z  <50-80  km  and  for  waves  of  A^izs)  ^10  km 
at  all  altitudes. 
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(a)  Saturation 


(b)  Infrared  Radiation  (IR) 


(c)  Molecular  Viscosity  (MV) 


Figure  16:  “Top-of-the-atmosphere”  contributions  to  gravity-wave  momentum  flux  dissipation  of  simulated  gravity 
waves  as  a  function  of  source-level  vertical  wavenumber  |m(Zi)|  and  horizontal  wavenumber  Ktot  at  z  =  200  km, 
where  close  to  100%  of  each  wave’s  momentum  flux  has  been  dissipated.  Contour  labels  are  percentages,  due  to  (a) 
breaking/saturation,  (b)  IR  radiative  damping  and  (c)  molecular  viscous  damping. 
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Figure  17:  (a)  Nq{z),  (b)  Noo(.z)  and  (c)  km{z)  derived  by  fitting  (A.l)  to  the  m)  estimates. 
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Appendix  A.  Least  Squares  Fits  to  ^(z,  m) 


We  have  fitted  the  seale-dependent  IR  radiative  damping  rates  in  Figure  4  to  various  functions 
to  facilitate  their  efficient  numerical  implementation  in  gravity-wave  models  and  parameteriza- 
tions.  Following  Zhu  (1993),  we  first  performed  a  nonlinear  fit  of  the  rates  at  each  height  z  to  the 
analytical  function 


r,j/z,  m)  =  Ao(z)  N^{z)  \  1 


tan  ^  [mikmiz)] 


(A.l) 


mlk^iz)  ) 

This  equation  fits  the  rates  as  the  sum  of  a  cooling-to-space  term,  No(z),  plus  a  layer  exchange 
term  based  on  the  analytical  solution  of  Spiegel  (1957)  for  an  infinite  homogeneous  atmosphere 
(Goody  and  Yung,  1989).  By  analogy  to  the  Spiegel  (1957)  solution,  km(z)  is  the  mean  absorption 
coefficient  and  No(z)  +  Nco(z)  is  the  damping  rate  in  the  transparent  limit.  Using  (A.l),  scale- 
dependent  radiative  damping  rates  for  the  reference  temperature  profile  and  constituent  mixing 
ratio  profiles  are  now  conveniently  parameterized  in  terms  of  a  lookup  table  of  reference  profiles 
of  No(z),  Ncc(z),  and  km{z),  which  are  plotted  in  Figure  17.  This  fit  has  the  advantage  that  it  is 
well-behaved  and  relaxes  to  physically  plausible  constant  values  at  m  =  values  outside  the 
=1-500  km  range  where  we  computed  raw  rates  from  our  Curtis  matrix  calculations. 

A  drawback  is  that  close  fits  to  the  original  data  do  not  occur  at  every  height  and  wavenumber 
using  this  function.  We  obtained  much  closer  fits  using  the  third-order  logarithmic  polynomial  fit 


r^j^(z,  m)  =  exp  [a(z)  -i-  ^(z)<A  +  c(z)<A^  -t  d{z)^^ , 
i/r  =  log  m  -  log  2n(500  km)“^ . 


(A.2) 

(A.3) 


While  this  fit  has  the  advantage  of  yielding  a  much  more  accurate  fit  at  all  heights  over  our  1- 
500  km  vertical  wavelength  range,  it  is  poorly  behaved  outside  of  this  wavelength  range  and  thus 
cannot  be  used  at  wavelengths  shorter  than  1  km  or  longer  than  500  km,  which  is  not  a  serious 
drawback  since  in  this  case  the  rate  at  the  limiting  boundary  wavenumber  is  imposed. 

Lookup  tables  containing  coefficients  of  both  fits  on  the  height  grid  Zk  are  provided  online  as 
supplementary  material,  in  the  form  of  two  text  files  designed  to  be  used  as  data  files  for  a  computer 
subroutine  that  reads  in  and  then  uses  these  coefficients  to  compute  damping  rates.  These  lookup 
tables  also  store  the  background  temperature  profile  T^izu),  so  that  fitted  rates  can  be  scaled  to  any 

local  temperature  profile  T{z)  using  (11)  and  (12). 
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Appendix  B.  Anelastic  Gravity  Wave  Model 


Our  one-dimensional  model  of  gravity-wave  evolution  with  height  is  based  on  the  anelastic 
gravity  wave  equations  of  Marks  and  Eckermann  (1995).  For  use  on  our  height  grid  zu,  we  set  the 
ray  time  step  At  to  Az/Cg^,  where  Az  =0.5  km  is  our  grid  resolution,  and 

-a>+m(l  - 

V  /  /T^  -1  \ 

—  T't  o  T’  v^'l) 

+  P  + 

is  the  wave’s  vertical  group  velocity.  We  next  define  a  complex  vertical  wavenumber  M  =  m  +  inii, 
such  that 


such  that  momentum  flux  varies  with  height  above  its  source  altitude  z  >  as 


F(z)  =  F(z0iS^(z)exp|-2  J  mi(z)dz  .  (B.3) 

S^(z)  is  a  scaling  factor  that  accounts  for  the  accumulated  loss  of  momentum  flux  due  to  wave 
breaking  (0  <  <S(z)  <  1),  and  is  described  below,  while  exp[-2  r  mi(z')dz']  is  the  time-integrated 

'^Zs 

loss  of  momentum  flux  along  the  ray  path  due  to  radiative  and  viscous  damping.  In  the  absence  of 
these  dissipative  processes,  both  terms  equal  unity,  whereupon  F(z)  =  F(Zi)  and  5F  jdz  =  0. 

We  initialize  momentum  flux  at  the  source  altitude  z^  by  specifying  the  wave’s  peak  vertical 
displacement  amplitude  ^(zj,  which  is  then  converted  into  a  corresponding  action  density, 

A  =  (B.4) 

4 

and  momentum  flux,  F,  via  (14),  where  the  polarization  factor 


P{z)  = 


(an  equivalent  derivation  in  terms  of  horizontal  velocity  amplitude  is  given  in  Appendix  B  of 
Marks  and  Eckermann,  1995).  Eqs.  (14),  (B.l),  (B.4)  and  (B.5)  imply  that  the  wave  amplitude 


varies  with  altitude  as 


^(z)  =  ^(z.) 


S'^{z)P{Zs)pbizd)CO'^{z)Cg^iZs) 

52(z,)E(z)pfo(z)m+(z.)Cg^(z) 
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i-r 


mi(z)dz  . 


As  the  wave  propagates  upwards,  at  each  height  Zk  we  check  whether  wave  amplitudes  exceed  the 
local  threshold  amplitude  for  wave  breaking  Isai- 


tiz)  >  lat 


a 

Iml’ 


(B.7) 


where  the  dimensionless  wave-breaking  amplitude  a  is  given  by  eqs.  (33)  and  (34)  of  Fritts  and 
Rastogi  (1985),  such  that  0  <  a  <  1  for  p  +  0  due  to  wave-induced  dynamical  (shear) 
instabilities,  and  a  — >  1  for  due  to  onset  of  wave-induced  convective  instabilities. 

We  adopt  the  standard  linear  saturation  hypothesis  used  in  most  parameterizations,  in  which  the 
wave  sheds  sufficient  momentum  flux  to  reduce  the  wave  amplitude  to  the  threshold  for  instability 
iX  =  tsat)^  which  is  achieved  by  setting 


.  ^sat(z) 

=  min  1, 

L  i(z)  J 

u 

(B.8) 

S(zk) 

K 

!=0 

(B.9) 
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